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Abstract 

Quantum Monte Carlo simulations of interacting electrons in solids often use 
Slater-Jastrow trial wave functions. The Jastrow function takes into account cor- 
relations between pairs of electrons. In simulations of solids, it is common to 
use a Jastrow function which is both homogeneous and isotropic. The use of 
a homogeneous and isotropic Jastrow factor is questionable in applications to 
systems with strong density inhomogeneities, such as surfaces and multilayers. 
By generalizing the original derivation of the RPA Jastrow factor for the ho- 
mogeneous electron gas ^ to inhomogeneous systems we derive an scheme for 
generating inhomogeneous and anisotropic Jastrow factors for use in nonuniform 
systems, from the non-interacting static structure factor and density of the sys- 
tem. We discuss aspects of our scheme and illustrate it with a first application 
to an inhomogeneous electron gas. 

1 Introduction 

In this paper we consider Slater-Jastrow trial wave functions, ^ = Dc"^ , where D is 
a Slater determinant and J the Jastrow factor that takes account of the electronic 
correlations. Such wave functions are employed in variational quantum Monte Carlo 
(VQMC) simulations of solids [0,0], where an approximate many-electron ground 
state is obtained by numerically optimizing an explicit parametrized trial wave func- 
tion. Remarkably better results are obtained from VQMC calculations using Slater- 
Jastrow functions than from Hartree-Fock (HF) or density-functional theory within 
the local density approximation (LDA). For example, cohesive energies obtained us- 
ing VQMC are typically an order of magnitude more accurate than those obtained 
employing HF or LDA [l-g. 

We use a Slater determinant of LDA orbitals and a Jastrow factor that includes 
pairwise correlation terms M(ri,rj). In simulation of solids u is usually taken to be 
both homogeneous and isotropic, i.e., u{ri, rj) is assumed to depend only on the inter- 
electronic distance = |ri — rj|. 

The long-range behaviour of u{ri, rj) is determined by the zero-point motion of the 



plasmons, both in homogeneous and inhomogeneous systems |10| while the short 



range behaviour is governed by the so-called cusp condition |jTT|. Gaskell minimized 



the variational energy of the uniform electron gas in the random-phase approximation 
and found a Jastrow function which combined both the short-range and the long- 
range behaviour. Ceperley used this Jastrow function in his variational Monte Carlo 



calculations of the spin-unpolarized homogeneous electron gas [12 



Here we generalize Gaskell's approach to inhomogeneous systems and derive an 
inhomogeneous and anisotropic Jastrow factor that depends on the non-interacting 
structure factor and density of the system. The Jastrow factor is represented as a 
double Fourier integral. The full interacting many-electron Hamiltonian H is expressed 
as a function of the plane-wave expansion coefficients of the Jastrow function and, 
within the random phase approximation (RPA) H is minimized with respect to these 
coefficients. This way the approximated ground state energy is obtained, and the trial 
wavef unction of the system can be constructed from the resulting inhomogeneous and 
anisotropic Jastrow function. In the homogeneous limit, our approach reduces to that 
of Gaskell. 

The rest of the paper is organized as follows. In Sec. 2.1. we describe the formalism 
used in Gaskell's original derivation of the Jastrow function for the homogeneous 
electron gas (HEG). Sec. 2.2 presents our generalization to systems with density 
inhomogeneities where, for the sake of simplicity, we consider systems with density 
variations in only one direction. Sec. 3 discusses the numerical results obtained for 
the IHEG, and Sec. 4 concludes. 



2 Formalism 



A Slater- Jastrow trial wave function is the product of a totally antisymmetric Slater 
determinant D and a totally symmetric Jastrow factor : \l/ = De^ . The orbitals 
used in D are obtained from LDA calculations. Here we consider spin unpolarized 
systems. The Slater determinant builds in exchange effects but neglect the electronic 
correlations caused by the Coulomb interactions. The most important correlation 
effects occur when pairs of electrons approach each other, and these may be included 
by choosing pairwise Jastrow factors of the form: 



where ti(rj, r^) 



J 



(1) 



u{Yj, Yi). So the trial wave function can be written as: 



^(R) = Dexp 



(2) 



Note that we include the i = j terms in the sum over i and j, which are omitted by 
many authors as they give a constant contribution for the HEG. (Unless otherwise is 
stated, we use atomic units throughout, i.e., = Ti = rrif, = 1.) 
The full interacting many-electron Hamiltonian is given by 



H = T + Vee + Ve 



exti 



(3) 



where T and V^e are the kinetic energy and the electron-electron interaction operators 
respectively and V^xt is an external potential. The expectation value of the energy 



with respect to the above trial wavefunction can be written as (the normahzed \1/ is 
taken): 



N 

E =< ^\H\^ >=Y.^^+ I drn{r)V,^t{r) - I rfm(r)K(r)+ < T2 > + < >, (4) 

1=1 



where Vs is the effective single-particle potential associated with the orbitals in the 
Slater determinant, are the energy eigenvalues which correspond to such one electron 
orbitals and n{r) is the electron density. < T2 > and < Vee > are defined as: 



< Ta >= J d i?^*(R)^(R) 



N 



Y: V^(E Uir„ r,)) ■ V„(E Uir„ r,)) 



m=l 



and 



< >= / rfRvi/*(R)vl/(R)(V^^). 



(5) 



(6) 



Here R = (ri, r2, . . . Tn) represents all electron coordinates. Several different types of 
Jastrow factors are commonly used, the best known are based on the RPA in the form 
due to Gaskell and the form based on the RPA of Bohm and Pines [§. Here we will 
briefly describe the Jastrow factor due to Gaskell for the homogeneous electron gas 
before presenting its generalization to inhomogeneous systems. 



2.1 The homogeneous RPA Jastrow factor 

For homogeneous systems the correlation term u depends only on the interelectronic 
distance rij and may be expressed as follows: 

M(ri,rj) = u{ri,) = Y.u{k)d^''^^ , (7) 

k 

where k is the modulus of the total wavevector. Also, the single-particle orbitals are 
a set of planewaves. Inserting the above expression in Eq.(p, the total energy in the 
random-phase approximation of T2 can be written as 0]: 

^ = E + 2iV' E kW{k)S{k) + 27^0 E ^^[S{k) - 1] (8) 

k<kp k k 

where NS{k) =< \E'|nknJ^|\E' > is the static structure factor of the uniform interacting 
electron gas and nk = are the number density operators, kp and Uq = N/V 

are the Fermi energy and the average electron densityi, respectively, being the 
total number of electrons in the system and V the volume. Gaskell assumed that the 
structure factor for the interacting electron gas is close to that of the non-interacting 
electron gas and, for a small Jastrow function ^(r) [| can be approximated by the 
perturbation formula 

^^^^ = i + Mk)So{ky 



"'^This corresponds to the high density hmit. 



where S'o(fc) is the structure factor of the non-interacting electron gas [0]. The varia- 
tional energy in the random-phase approximation (Eq. ^ can then be minimized with 
respect to u{k) to yield: 



u{k) 



-1 1 



1 IGvrnn 
+ 



1/2 



2.2 The inhomogeneous RPA Jastrow factor 



(10) 



Here we consider the case of an inhomogenous system translationally invariant in the 
plane perpendicular to the z axis. In this case u can be written as: 



u 



where r = (ry, z) and k = (q, k^). u^^ are the two-dimensional Fourier coefficients 
of the correlation function u: 



Uk 



(12) 



We use the following definition of the number density operator: nq,/;^ = Y.i e*'^ '"ii»e 
so that the Slater- Jastrow trial wave function is given by: 

\E' = Dexp 



(13) 



2.2.1 Deriving an approximation for the inhomogeneous structure factor 

In this section we develop an approximation based on the RPA to obtain the interacting 
structure factor in terms of Fourier coefficients of the above Jastrow factor. Let us 
define: 



^(A)=L>exp -A^ ^ Uk^„k,,{(l^n^MKk. 
The expectation value we need to calculate is of the form: 



(14) 



/(A) 



<^(A)Kfc.,n*, |^(A)> 



< ^(A)I^(A) > 
with A = 1. With the above trial function this becomes 



J(A) 





.,exp 


-2AEq2 T.k,^,k,^ Uk,^,k,^{q,2)nci2,k,^ 


"'q2fcz4. 


J DD*exp 


-2A Eq2 Efc,3,fc,4 «fc.3,fc.4(q2)^q2,fc.3 5^q2fc,^ 





(15) 



(16) 



Following Gaskell we use the identity: 

1 1 



rX Q 


1 


lo d\' 





dX' 



(17) 



to obtain a perturbative approximation for /(I). Taking the derivative inside the 
integral with respect to A and expanding in a Taylor series the leading term is 



<0|("qi,fc.l";ifc,,)'|0>-<0|"qi,fczi"q,fc,^|0> 



Mfcz3,fcz4(q2) 

<0|"q2,fcz3";2fe,,"qi,fezi'^;ifc,JO> <0| "q2 ■ ^zg fe,, 1 0><0| {«qi ,fezi fc, J 0> ) 



3 q2'°Z4 qi''^zi qi/Cz2' ' q2."-Z3 q2fcz4 ' ' ' " ML-^Z]^ qi fcz2 ' 

qi'=z2' ' qi.Kzj qifcz2' 



(18) 

where < 0|A|0 > denotes the expectation value of an operator A with respect to 
\l/(0) = D. Following Gaskell, we use the following approximation: 



< 0\{h^,k.n*^,J\0 >= 2 [< 0\{h^,k.n*^,J\0 >]' + 0{1/N), (19) 

and ignore correlations between different h's in the second term. The resulting ap- 
proximation for /(A) becomes 

'^^^ = l + A/(OHl^,..^(q)- ^2°) 

The corresponding approximation for the Fourier coeffecients of the structure factor 
is obtained form the above equation for A = 1, yielding 



Sk.„k.M) = , , '"tAoo — (21) 



2.2.2 Expectation value of the energy 

In order to obtain the expression for the total energy of an inhomogeneous system we 
use the RPA approximation 

< '^qi,fczi^q2,fcz2^qi+q2,fcz3 > = < '^qi,fczi'^q2,fcz2 '^qi+q2,fcz3 > • (22) 

In the case of a system with density variation along z we have: 

^e, ^- / if qi + q2 = , . 

< ^qi+q2,fcz >- I Otherwise ^^^^ 

where p^z represent the Fourier coefficients of the electronic density. 

Making these assumptions and employing the expression for the structure factor 
obtained in the previous section, the terms < T2 > and < V^e > which appear in the 
expression for the total energy Eq.(^ can be written as: 

< T2 >= 2VA^Eqi,fc,^,fc,2,fc,3,fc,4Mfczi,fez2(9l)Wfez3,fcz4(9l) [(?l)^ " ^^l^^s] 

XPfczi+fcz3 (^.Z2,-fcz4(?l)) (24) 

< Vee >= noEci,k. [('^^•--'^^(^)) - 1] • 



The above equation combined with our approximation (Eq.^Tj) for the structure factor 
leads to an expression for the total energy in terms of the Fourier coefficients of the u 
function, the electron density and the non-interacting structure factor. We take both 
the density and the non-interacting structure factor from LDA calculations and mini- 
mize E with respect to the Uk^,k^,{l) coefficients, which then can be used to construct 
the trial wave functions for VQMC calculations. 



3 Results 

This section presents the numerical results obtained for uniform and non-uniform spin- 
unpolarized A^— electron systems (A^ = 38) contained in a simple cubic (SC) simulation 
cell and satisfying periodic boundary conditions within the cell. Density modulations 
in z are induced by subjecting the system to and external potential of the from 

Vext{r) = Vocos{qz) (25) 

with Vq = 0.1 a.u.. q = O.Gkp is a reciprocal lattice vector of the SC simulation cell. 
The average density of the system corresponds to r<j = 2. The electron density is 
uniform in the x — y plane and varies in th z direction. 

The inhomogeneous u(ri,r2) term is represented by fixing the position ri of the 
first electron, while moving the second electron, r2, along the z direction on a line 
passing through ri. 

3.1 Homogeneous system 

In Eq.(|^) we saw that u{zi, Z2) may be expanded as a Fourier series. In Fig. |I| the 
convergence of the energy with respect to the cutoff wavevector, kcutoff, used in the 
Fourier expansion of the u factor is shown. The typical cutoff needed is small enough 
to be computationally feasible, so that such a representation for u{zi, Z2) is found to 
be very useful. 



kmiiir (in units of kp) 



Figure 1: The convergence of the total energy as a function of the cutoff kcutoff used in the 
truncated Fourier series representation of the u function. The results are for = 38 and 
Tg = 2. kcutoff is given in units oi kp. 

Looking at the expression for the energy we see that for u{k) = 0, < T2 >= 
and < Vee > is just the exchange energy of a HEG. For non-zero u{k) factors 
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ec{Gaskell) 
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-0.04677 


-0.058 


3 


-0.04175 


-0.05 


4 


-0.03809 


-0.0445 


5 


-0.0352 


-0.0405 



Table 1: Correlation energy per electron in a.u. 

< T2 > -|- < Vee > includes both exchange and correlation energies. The correlation 
energy is defined as: 

= [< T2 > + < Vee >] - [< T2 > + < Vee >]«(fc)=0 (26) 

Table 1 presents the values for the correlation energy obtained in this work in com- 
parison with those obtained by Gaskell|I| . Difi^erences due to finite size effects have 
to be taken into account, as we have considered a very small system with only = 38 
electrons. 

The results for the energy obtained with the uniform system previously described 
and using the homogeneous RPA u factor given in Section. 2.1. are presented in Table 
2. 

3.2 Inhomogeneous systems 

The inhomogeneous RPA jastrow factor considered here is the one derived in Section. 
2.2. The energies derived are presented in Table 2. Exchange-correlation energy is 
defined as: 

E.e =< T2 > + < Vee > 'Eh, (27) 

where represents Hartree energy. 





kin 






^tot 


I) 


0.2807 


-0.157 


-0.04677 


0.1235 


11) 


0.2658 


-0.1193 


-0.0061 


0.144 



Table 2: Non-interacting kinetic energy e°j^ = j^Y^k^"^ f^i exchange-correlation energy 
€xc = Exc/N, correlation energy Cc = Ec/N and total energy etot = E/N per electron in a.u. 
for I) an uniform system using the homogeneous RPA u factor, II) a weakly inhomogeneous 
system using the inhomogeneous RPA u factor. The cutoff wavevector used was equal to 
2.5/cF. 

In Fig.^ we plot the inhomogeneous RPA u term, for the three different positions 
of Fi indicated in Fig. |]b which presents the LDA density and the average elec- 
tron density. The homogeneous Jastrow factor is shown as well, which depends only 
on the relative distance between the electrons and not their individual coordinates. 
Anisotropy of the Jastrow factor can be observed. This is such that the u function is 



3 6 9 -50 -25 0.0 2.5 5. 

z (atomic units) z (a.u.) 

Figure 2: a) LDA electronic density of the system. The horizontal line represents the 
average electronic density of the system, b) Inhomogeneous RPA u term for the three 
different positions of the fixed electrons shown in a), together with the homogeneous result 
for the same system. 

stronger on the side where the density is lower. This can be explained by the fact that 
the RPA screening is more effective where the electron density is high, so that the u 



function is weaker on the high density side |]T0| 



4 Conclusions 

We generalize the RPA treatment of the homogeneous electron gas by Gaskell to the 
inhomogeneous case. As a result we obtain a Slater- Jastrow trial wave functions con- 
taining an anisotropic and inhomogeneous Jastrow factor. The latter may be obtained 
from the non-interacting structure factor and density of the system. 

Results for the energy and the inhomogeneous correlation u factor using an elec- 
tronic system subject to a sinusoidal external potential are presented. We compare 
these results with those obtained using the homogeneous u factor by Gaskell, and 
we clearly observe the anisotropy introduced by the inhomogeneous u factor. Future 
work will focus on testing the performance of our inhomogeneous Jastrow factor using 
variational Monte Carlo and in further refining our approximation of the interacting 
structure factor. 
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